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One-dimensional systems have unusual properties such as fractionalization of degrees of freedom. 
Possible extensions to higher dimensional systems have been considered in the literature. In this 
work we construct a mean field theory of the Hubbard model taking into account a separation of 
the degrees of freedom inspired by the one-dimensional case and study the finite-temperature phase 
diagram for the Hubbard chain and square lattice. The mean field variables are defined along the 
links of the underlying lattice. We obtain the spectral function and identify the regions of higher 
spectral weight with the fractionalized fermionic (spin) and bosonic (charge) excitations. 



I. INTRODUCTION 

The Hubbard model is one of the simplest models ac- 
counting for interactions between electrons and has been 
used to describe various strongly correlated materials. It 
is parameterized by only two constants, a kinetic term 
scaled in the tight-binding approximation by the hop- 
ping amplitude t and an on-site interaction U modeling 
the Coulomb repulsion felt by two electrons of opposite 
spins occupying the same site. 

The exact solution of the one-dimensional Hubbard 
model by the Bethe ansatz [T| involves composite degrees 
of freedom that correspond to different rapidity branches. 
In addition to a cO charge-momentum rapidity, there are 
sets of (oc,u) rapidities [2]. The general rapidity branch 
label av is such that a — c, s and v = 0, 1, 2, ... for a = c 
and v = 1,2,... for a = s. The (c, v) and (s, v) rapidi- 
ties are associated with the charge and spin degrees of 
freedom, respectively. For electronic densities n < 1, the 
ground state has finite occupancies for the charge (c, 0) 
and spin (s, 1) branches only. The relation between the 
original electrons and the entities that describe the eigen- 
states of the Hubbard model is however complex and only 
recently some light has emerged [3]. 

Electron double occupancy is a good quantum num- 
ber (it is conserved) in the limit of U — » oo but for 
finite values of U/t it is not conserved. However, it is 
possible to define new fermionic operators, associated 
with fermionic objects called rotated electrons, through 
a canonical transformation, V, such that the double oc- 
cupancy of these rotated electrons is a good quantum 
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number for all finite values of U/t [3— 6J . In terms of 
the rotated electrons, it is a consistent interpretation of 
the Bethe ansatz states to describe the various branches 
(rapidities) in terms of, first, a separation into empty, 
double-occupied and singly occupied sites. The empty 
and double-occupied sites are called 77-spin 1/2 holons [3 
and the singly-occupied give rise to a charge part that 
originates the (c, 0) particles and to a spin part which 
originates the spin-1/2 spinons. The fy-spin projection 
1/2 (and —1/2) holons correspond to rotated-electron un- 
occupied sites (and doubly-occupied sites). The spinons 
of spin projection ±1/2 refer to the spins of the rotated 
electrons which singly occupied sites. Second, these are 
paired in singlets (s, 1) or in pairs of pairs of singlets and 
so on (s, v). The holons are also paired in such a way that 
empty sites and doubly-occupied sites are paired (c, v) . 
These results apply to the rotated electrons and not to 
the original electrons. They are however related by the 
above mentioned canonical transformation. At high val- 
ues of U they are very close and identical when U — > 00. 
For many practical situations, such as the calculation of 
various correlation functions, it is a reasonably good ap- 
proximation to consider the rotated electrons as similar 
to the original electrons [THTT]. 

There are transformations in the literature that pro- 
pose a similar decoupling of the electronic degrees of free- 
dom. This can be seen for instance in j 12 or in |13j . 
The main motivation was the study of either the large- 
U limit in the Hubbard or Anderson models |14j with 
the intent to control in an efficient way the projection 
to states where double occupancy is restricted (as in the 
t — J model) but considering a finite value of U instead of 
the extreme case of infinite U, usually taken care of by a 
single slave boson [15]. In the Zou- Anderson transforma- 
tion (ZA) each physical electron is mapped into the one 
particle sector of a set of four particles, two fermions and 
two bosons. The two bosons may be chosen as spinless 
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particles and represent the empty and doubly occupied 
states of each lattice site; excitations of these two de- 
grees of freedom are called holons and doublons and can 
be interpreted as carrying no charge or — 2e respectively. 
The fermions then represent singly occupied states with 
spin ±1/2, excitations on this sector are called spinons 
and carry charge — e. However, the charges may be de- 
fined differently, as considered in the original paper by 
Zou and Anderson [T3]. There is also an exact transfor- 
mation that explicitly includes spin and charge separa- 
tion [TB] introducing two sets of operators as quasicharge 
(fermionic) and quasispin (spin-like). The representation 
of the quasispin operators leads however in general to 
bilinear terms in bosons or fermions and therefore this 
leads to six operator terms in the Hamiltonian. 

The ZA mapping reverses the role of the interacting 
and kinetic terms in the Hamiltonian. The interacting 
Hubbard term becomes quadratic in the ZA particles and 
the kinetic one is transformed into an interacting quartic 
term that couples particles along the lattice links. This is 
particularly useful to study the strong interacting (large 
U) regime where the kinetic term is treated as a per- 
turbation. The price of this transformation is the ap- 
pearance of an on-site constraint which assures exactly 
one particle per lattice site. In the mean field (MF) ap- 
proach this translates to an on-site Lagrange multiplier. 
Generically slave-particle methods induce new unphysi- 
cal symmetries of the Hamiltonian written in terms of 
slave particles, that are called in this context gauge sym- 
metries. In this particular case the symmetry group is 
U(l). 

We note that the representation introduced by Zou and 
Anderson has been used to explicitly obtain an exact 
solution of the Hubbard chain in the large U limit in 
a much simpler way as compared to the Bethe ansatz 
|17j . Also it has been used to study the stiffness of the 
one-dimensional Hubbard model in a way equivalent and 
complementar to the Bethe ansatz solution |18| . 

In the large U limit it is very costly to create doubly 
occupied sites. It is usual to consider d = jTS] (projec- 
tion to subspace of no-double occupancy) or to perform a 
canonical transformation to eliminate transitions to dou- 
bly occupied sites [13| . This type of procedure leads to 
a treatment similar to the t — J model, appropriate near 
half filling in the vicinity of magnetic order for a square 
lattice. 

It has been argued recently that one should instead 
integrate out the high-energy scale to obtain an effec- 
tive low-energy theory which has been shown to contain 
a charge 2e bosonic mode [20 22 . This mode may be 
bound to a hole, providing a possible explanation of the 
pseudogap observed in high-T c materials. In this work 
we will maintain the full structure of the transformation 



between the electron operators and the auxiliary parti- 
cles introduced in the Zou- Anderson transformation. We 
will be focusing on the finite energy (finite temperature) 
properties where we will find and study phases with non- 
standard correlation functions associated with the frac- 
tionalized degrees of freedom. We will obtain the phase 
diagrams for one- and two-dimensional systems and cal- 
culate various correlation functions as well as the spectral 
function. 

We will use nonlocal decouplings of the auxiliary op- 
erators leading to link variables that can be associated 
with nearest-neighbor spin singlets or bound-states be- 
tween empty and doubly occupied sites in a way close 
to the Bethe ansatz solution and also suggested by the 
treatment of Ref. imposing at MF level the existence 
of these states. A simplified treatment of link variables 
was introduced in Ref. [23 and the significance of short- 
range correlations between empty and doubly occupied 
sites was, for instance, determined in Ref. [53]. We note 
that bond variables appear naturally in extensions of the 
Hubbard model too, for instance nearest-neighbor inter- 
actions, or bond correlated hoppings [25j. Interestingly, 
many results have been obtained for systems where the 
hopping between singly occupied sites and empty and 
doubly occupied sites is eliminated, implying that dou- 
ble occupancy is a good quantum number |26H28| . 

The paper is organized as follows: in section II we dis- 
cuss the methods used to construct the MF solution. In 
section III we discuss the MF phase diagram for the Hub- 
bard chain and the square lattice and in section IV we 
discuss the spectral function for both cases, interpreted 
in terms of the fractional excitations here considered. Fi- 
nally in section V we discuss the results obtained and in 
the appendix we present a list of MF solutions for the 
square lattice. 



II. METHODS 

We start from the U(l) representation of the Hub- 
bard model introduced by Zou and Anderson (ZA) |13| . 
<4 a — e r sj, a + ad\ s r - a which maps a spin-1/2 fermion 
(the physical electron) to the one-particle subspace of 
four fields e r ,d r ,s r>a —±x- In this work bosonic e r and 
d r fields (where r labels the lattice sites) are considered 
corresponding to the annihilation of empty and doubly 
occupied sites, s r CT= ±i are fermions that carry the spin 
degree of freedom. The opposite choice of statistics (e, d 
fermions and s±i bosons) is also possible leaving the 
mapping unchanged. The Hamiltonian of the Hubbard 
model is given by 
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H - fin T = -t ^2 cl a c r+ s,a + cl +S a c r ^ + U^ j n r ,in r ,-i - \i ^ n r . a , 

r r,a 



(1) 



r',5>0,CT 



r,(5>0 



+ (d r e r+ 5 + e r d r+ s)^ {s r ,lS r +8-l — Sr,-lS r +<5,l) + ( s r,lSr+5,-l ~ Sr,-lSr+5,l) (d r e r +5 + e r d r+< 5) 



(2) 



where 5 is a directed lattice vector connecting nearest 
neighbor sites, /i is the chemical potential and Ut the 
total number of electrons. 

Using the ZA mapping the partition function is 
given in a path integral formulation by 
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r,S r 



where 



B, 



1,3 



d]<k-e\ 




d\e\- 




eidj + ejC 


k 






a i-i*j,-i 




4-i4i - 4i s i,-i 
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4,-i s «.-i + s i,i s «,i 



(3) 



are respectively bosonic and fermionic matrices, A r is 
an on-site real field inserted in order to project to the 
one-particle sector n^ A = d\.d r + eJ r e r + Y^, a 4 o- s r,<r = 
1 with/ DX = f^ UrUT i^, and V is the 
total number of lattice sites. Using the identity: 
fd f i(Q\Q)e- Tr W-Q+ c> -Q + Q' , - A } = e T ^- A \ with 
dQ i4 dQ M _ and 



Q' 



Qi 



[Q]i j — Qi.j(Qj,i — Qi j)i we introduce a 2 x 2 matricial 
Hubbard-Stratonovich (HS) field in order to decouple the 



fermionic and bosonic terms: 

Z = J DXDeDdDs a DQ e~ ^ drh (6) 
L = 4(<9r + X r + U - n)d r + Y e t(dr + X r + fi)e r 

r r 

+ Y S i,<r(dr + X r )s r<(7 - Y "V 
r r 



r.S 



+ ^v^Tr [Ql s .F r ^ 

r-\-5 ^r/r+^-Qr.i 



(7) 



The additional gauge freedom introduced when writing 
the Hubbard model in this particular slave-particle form 
is U(l) and is implemented by the operator U((j)) = 
Y[ r e~^ r ' lr ''. Such gauge transformation changes the ZA 
particle fields by a site dependent phase a r —> e" l< ^ r a r 
(a = s±,d, e) which gives the simple transformation rule 



for the matrices B 



r tr+ s -> e~'"- Vr + s .B,, ir+5 .e 4<T ^'' and 
In order to leave the 



Lagrangian invariant, this transformation also induces a 
gauge transformation in the HS fields 



Q 



r,6 



Q, 5 .e^+*, 



(8) 



defining also 6 on-site gauge invariant quantities: 

\9r, S ]i tj [Ql,i\ tj aild [QT,s] id [QrAj^^ 2 )- 

A MF treatment of this Lagrangian can be justified 
introducing N copies of the ZA particles in order to per- 
form a 1/N expansion that coincides with the usual MF 
approximation at zero order and organizes the follow- 
ing corrections in numbers of loops fT^]. However, as in 
SU(N) generalizations of spin-1/2 models, this param- 
eter is rather unphysical and this approach will not be 
explicitly pursued here. The MF approximation is ob- 
tained varying the free energy with respect to the Q and 
A fields: 



QrA T ) 

i 



ts (F r , ff (r)) 
(n ZA (r))) , 



(9) 
(10) 
(11) 
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where ( ) stands for the MF average. In order to study 
this set of equations a time and space translational in- 
variant ansatz Q r g(r) = Q s was imposed. The saddle- 
point values of the HS fields describe hopping and pairing 
terms: 

XF,r,8 = ( S Y-{-8 t l s r,l ~t~ S r+8. — l s r, — 1 

^F,r,S = (s r +S,lS r ,-l — S r +S,-lSr,l) l 

J 



XB,r,S 



d r +S e r+S C r 



&B,r,6 — (d r +se r + e r+ sd r ) Q , 
and Eqs. ( 9|[T0 I are equivalent to 



(12) 



( v f A f 

r.S - ~WH . ' 

\ &F,rJ XF,r,S 



Q 1 



r,6 



, , XB,r,S &B, h S 
S I As,r,5 X B .r,S 



(13) 



Note that the MF values obtained for and Q are 
not complex conjugated of each other; this is crucial 
in order to interpret the zero order Lagrangian as com- 
ing from a hermitian MF Hamiltonian, in the extended 
Hilbert space, as noticed in |2"31 129]. This corresponds 
to the analytic continuation of the Q fields and also oc- 
curs for A. Care should be taken, however, when consid- 
ering fluctuations of these fields around their MF val- 
ues: as the A fluctuations are purely imaginary even 
if the MF value is real, the conjugated fluctuations of 

SxF.r,s are $X B rS an< ^ no ^ ^XprS as t ne MF treatment 
could suggest. The filling fraction of the electrons is im- 
posed as usual requiring the chemical potential to satisfy 
x = 1 — (ut) = 1 — \7/3^m m where x is the hole doping. 

Assuming translational invariance of the MF solutions, 
the MF Hamiltonian can be brought to a diagonal form 

in k space. Defining b k = jdfc, e-k,i, d\. <v e^_ k and 

Sk = |sfe,i) s-k,-i, slfe^xj the Bogoliubov trans- 
formation diagonalizing the Lagrangian is given by 



a k = R b (k).b k 
f k = R s (k).s k 



(14) 



where f a =± and cii =e .d are the fermionic and bosonic Bo- 
goliubov transformed fields . In these new variables the 
Hamiltonian is given by 



H 



a k,i e i,k a k,i + ©B,fc 

k,i=d,e 

+ fk,* £ F,kfk,* + ®F,k 



s 



where the single particle energies e and the energy shifts 
are given by 



£F,k = ^ \^>F,k - u F -k + \, J-'V.a'V./. + U'r. /, ^j.-.i, )' 

©F,fc = — 2 \- £p : k + £ F-k 



+ XB(k) - x B ( k ) + XB(-k) - x\{-k) 
1 



,t 



-•I.I, : 2 I Ud ' k ~ U e-k + V (W e ,_fe + Wrf.fe) - 4:5 B ,kSB,k 



-.J. -~ 7i ( We,k — Vd,-k + V ' {^e,k + ^d-kf - ^B,k&B,k 



©S,fc — 2 \- ed ,k + e e~k 



with 



+ XF{k) - xU k ) - XF{~k) + xU-k) 



u>d,k = U-n+ X + XF(k) + x F ( k ) 

w e ,-* = A* + ^ - XF(-k) - x F (- k ) 

6 B , k = A F (-fc) + A F (fc) 

u F ,k = A + xs(fc) + X B (k) 

S F ,k = A B (-k) + A B (k) 



(16) 



(17) 



and A(k) = £ 5=out tse tk5 A s (A = X ,A). Explicitly the 
Bogoliubov rotation matrices are given by 



R, = 



Rh = 



( u F , k v F ,k \ 

u Fk -v Fj k 

v Fk u Fk 

\ -v Fk UiT, fc / 

( ilfii V B ,k \ 

its,fc w Bife 
v B ,k u Btk 
\ v B , k j 



(18) 



with 

uf/bm 

VF/B.k 

r F (k) 
r B (k) 



V 1 + T F / B {k)T F / B {k) 

tf/b 

V 1 + T F/B(k)f F/B (k) 



(uj Fk + u>F,-k) - V (uF,k + u F -k) + 4:S Fk 5 Fk 

2$F,k 



(w d ,k + w e ,_/t) - V (wd.fc + w e ,-fc) - 4:5 Bik S B ^ k 
2S B h 



The above treatment should be valid for arbitrary val- 
ues of the interaction parameter U once the double occu- 
pancy described by the d bosonic field is fully taken into 
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Figure 1: Mean Field phase diagram in the x — T plane as a function of U for a chain with 200 sites. The colors represent the 
different types of solution minimizing the free energy: (1) \ 7^ 0, A = (Red); (2) x = 0, A / (Orange); (3) x 7= 0, A / 
(White); (4) Xi A = (Black). The calculated points are placed in the nodes of the finite mesh of the values of x and T 
considered, the colors plotted between nodes are interpolated. 



x=0.1 




netic spin-spin correlation functions for non-frustrated 
lattices. In this work no boson condensation is consid- 
ered away from T = 0; this corresponds to a fully ID or 
2D model of holons e and doublons d UHl . 



Figure 2: Mean Field phase diagram in the U — T plane for 
x — 0.1 for a chain with 200 sites. The colors represent the 
different types of solution minimizing the free energy as in 
Fig. [1] 



account. This should be of great importance near half- 
filling because in this regime (c^ii) is not small compared 
to (ete). 

Contrary to the t — J model no canonical transforma- 
tion was performed at this stage in the physical electrons. 
Note that discarding the d field at this stage would yield 
diagonal F and B matrices and the subsequent MF de- 
coupling would miss the phases of the ZA fermions with 
non-zero pairings A. This kind of phases could however 
be considered if one adopts a variational procedure (Ref. 
|30|). The t — J model also presents an S • S term not 
present in Eq. ^ which is responsible for sub-lattice 
magnetization in the ground state of the square lattice 
Hubbard model near half-filing. Anti-ferromagnetic cor- 
relations will nevertheless be generated if one integrates 
out the d fields before doing the MF decoupling; however, 
such procedure leads to six-body coupling terms and will 
not be considered here. Even if Eq. ^ can not pro- 
duce sub-lattice magnetization we expect antiferromag- 



Eqs. (9fTT I together with the fixed doping condition 
were solved numerically in one and two dimensions for 
several values of U/t (typically in the range of 2 to 6) , 
for different values of doping x ~ 0, ...,0.6 and temper- 
ature. For the 2D case no particular symmetry based 
ansatz is implemented leading to some peculiar phases. 
Before discussing some of the MF solutions in detail a 
few remarks are in order: 

(i) Some care was taken diagonalizing the fermionic 
and bosonic quadratic Hamiltonians since the hoppings 
and pairings found were in general complex numbers. 
However some of the converged solutions lead to non- 
Hermitian Hamiltonians. These solutions were discarded 
as "unphysical". However, for the slave-particle approach 
the only physical constraints are for the composite elec- 
tron operators and so maybe some of these solutions 
could be physically meaningful. 

(ii) The convergence of the solutions was quite difficult 
for some regions, specially near zero doping and for very 
low temperatures. For two dimensions a few points in 
the U — x — T space were tested against finite size effects 
running the calculations in a 32 x 32 lattice (instead of 
16 x 16) with only small quantitative changes in the re- 
sults. However near regions of phase competition small 
finite size effects can change the global minimum yielding 
qualitatively different results. 

(iii) For real physical systems described approximately by 
the Hubbard Hamiltonian other small coupling terms are 
expected that can qualitatively change the phase diagram 
locally. Different stable solutions that are not global min- 
ima will be considered elsewhere and are listed in the 
Appendix. 
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Figure 3: Left panel: Free energy for ID, U — 4, x = 0.1 as a function of the temperature. The Black, White and Red dots 
follow the same color code as the phase diagrams. Central panel: Fermionic and bosonic order parameters as a function of 
temperature for each phase separately. Right Panel: Minimum energy as a function of temperature for the fermionic and 
bosonic bands. 



III. MEAN FIELD PHASE DIAGRAM FOR 
HUBBARD CHAIN AND SQUARE LATTICE 

A MF treatment of the original electron problem is 
not expected to be a good approximation for the one- 
dimensional case. However, it is interesting to construct 
a MF treatment in terms of new operators and to test 
the differences with the existing exact results and the 
domain of validity of the present treatment. It also per- 
mits to have a generic idea of the phase diagram and of 
the different phases present in higher dimensions. 

In the case of the Id Hubbard model we present some 
of the MF solutions obtained numerically solving the MF 
Eqs. (9fiTl on a 200 site chain. The results were ob- 
tained as follows: for a chosen point of the parameter 
space U — T — x we generated several random trial solu- 
tions and used them as a starting point to our numerical 
routine. We obtained several different solutions unrelated 
by gauge transformations; however, just 4 of these solu- 
tions are relevant for the range of parameters considered 
here. The other solutions have very high values of the 
free energy or are "unphysical". These 4 solutions were 
extended to the rest of the U — x — T parameter space 
using as initial conditions a nearby converged solution. 
The four physical MF solutions differ by the existence of 
non-zero hopping and pairing terms. Fig. [I] shows the 
phase diagram x — T for different values of U and Fig. 
[2] the phase diagram U — T for x = 0.1. In these figures 
the results are presented on a finite mesh of points in the 
x — T or U — T planes. 

For the two dimensional case we present the most sta- 
ble MF solutions obtained numerically by solving the MF 
Eqs. [9][TT| on a 16 x 16 square lattice. The results were 
obtained as follows: for a chosen point of the param- 
eter space U — T — x we generated 1000 random ini- 
tial conditions and used them as a starting point to our 
numerical routine. From this 1000 initial conditions 30 
different solutions, non gauge equivalent, were found to 



converge. After this first step the 30 different solutions 
were extended to the rest of the U — x — T parameter 
space using as initial conditions a nearby converged so- 
lution. This procedure is tedious since the convergence 
of the solutions is not always easy and some of them do 
not converge even if the difference in the parameters is 
small. Finally for each point in the U~x~T space (again 
defined on a finite mesh) the solution with smallest free 
energy was found. Note that we did not impose any par- 
ticular symmetry to solve the MF equations, the only 
assumption being translational invariance (Q r s = Qg), 
in order to diagonalize the system in momentum space. 
That fact explains the proliferation of solutions of the 
MF equations. 



A. Description of the Phase Diagram 

Generically the phases found solving the MF solutions 
are characterized as follows: 

- Phase (1) Conducting phase characterized by x ^ 
0, A = (Red): the spins are gapless and the charge 
degrees of freedom present a gap of the order of the tem- 
perature, which closes at T = 0. 

- Phase (2) x — 0) A ^ phase (Orange), gapped 
for both degrees of freedom. Since it appears near x — 
it is tempting to identify this phase as an insulating 
antiferromagnet . 

- Phase (3) x 7^ 0, A ^ phase (White): precursor of 
the superconductor, in this phase there exists spin singlet 
formation but the charge motion is incoherent since no 
condensation was allowed. If we had imposed (efc=o) = Z 
this phase would split in two sub-phases analog to the 
pseudogap and superconducting phases in |19| . 

- Phase (4) High energy phase (Black) : is an incoherent 
phase where all correlations are zero. 
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B. Characterization of the Phases 

Our numerical results in one and two dimensions show 
that the x — 0, A 7^ (Orange) phase is dominant 
in the low doping region up to a temperature that 
decreases with U. In the ID case this phase extends 
to the under-doped region and interfaces with the 
\ 7^ 0, A = (Red) phase, present at higher doping, 
by a small finite-energy region where ^ / 0, A 7^ 
(White phase). In 2D the x = 0, A ^ (Orange) phase 
is numerically unstable and we could only find it for 
zero doping. The x / 0, A 7^ phase appears also 
in the high doping regime at low T. In ID the size of 
this region grows clearly with U, however, in 2D this 
is not clear but is definitely present at low energy. At 
very low T the most stable solution is the Red phase 
(X ^ 0; A = 0) except at half-filling. 

One-dimensional case 

Figure [3] shows the free energy of the Hubbard chain 
for the four MF solutions when U — 4,x — 0.1 asa 
function of the temperature. At low temperature the 
X 7^ (Red) phase has a lower free energy and there is a 
first order phase transition with increasing temperature 
as the x / 0;A / phase becomes less energetic. At 
higher temperatures two second order phase transitions 
occur: first, for T ~ 0.23, x decreases to zero as the 
X 7^ 0; A 7^ joins the x = 0; A 7^ solutions; the second 
phase transition occurs for T ~ 0.27 when A vanishes 
(see Fig. |3j central panel). 

In the right panel of the same figure we show the mini- 
mum energy of the various bands as a function of temper- 
ature for the various phases. At low temperature where 
the Red phase is the most stable the fermionic bands are 
gapless. This phase is always gapless up to the point 
where it merges with the fully incoherent (black) phase. 
In the other phases the fermions are gapped (note that 
the spin spectrum in the Orange and White phases is al- 
most the same). On the other hand, the bosonic (charge) 
spectrum is always gapped except at zero temperature. 

Considering for instance U — 4, T — 0.01 we may as 
well analyze the results as a function of doping. At zero 
doping both x / 0; A ^ and x — 0; A 7^ coincide, 
this degeneracy is lifted for finite doping (in a very narrow 
region) and the x 7^ : A 7^ phase presents the minimal 
value of the free energy. From x~0.1toa;~0.5 there is 
an intermediate phase x 7^ 0; A = delimited by two first 
order transitions. For higher doping the x / 0:A ^ 
regains the minimal value of the free energy. 

The hopping and pairing correlation functions 

XF{r) = ^4,l s 0,l + 4,-l S ) 
AfW = (flr,l«0,-l — Sr,-lS0,l) o , 

Xb{t) = (4 do - e],e Q ) Q , 

Afl(r) = (d r e + e r do} , (19) 

are shown in Fig. [4]at two doping values x = 0.1, x = 0.6. 




2 4 6 8 10 5 10 15 20 

r r 

Figure 4: Hopping and pairing correlations as a function of 
the distance in ID, computed for the three non-trivial MF 
phases with U = 4. Red(full): Phase \ 7^ 0, A = computed 
for T = 0.10,2; = 0.1. Red(dashed): Phase % + 0, A = 
computed for T = 0.15, a; = 0.6. White: Phase x / 0, A 
computed for T = 0.19, a; = 0.1. White (dashed) : Phase 
X + 0, A computed for T = 0.08, x = 0.6. Orange: Phase 
X = 0, A / computed for T = 0.24, x = 0.1. 



Even if these correlation functions are not gauge invari- 
ant they can be quite useful to characterize the differ- 
ent phases. In particular one can clearly see the differ- 
ence between the two disjoint x 7^ 0; A 7^ (White) ((3) 
and (3')) regions. At lower doping 1 = 0.1 we consider 
the three non-trivial solutions as a function of increasing 
temperature and at the higher doping (x = 0.6) we con- 
sider the White and Red solutions. Both for the fermion 
and the boson hopping correlation functions the corre- 
lation length increases as the doping increases. Partic- 
ularly the bosonic correlation function has a large cor- 
relation length. Analyzing the correlation length of A^ 
one clearly sees a long range correlation in the high dop- 
ing regime possibly precursor of Bose-condensation and 
superconductivity. In the low doping region both the 
bosonic and the fermionic correlation functions have a 
smaller range consistent with a spin gapped state. In this 
regime the two correlation functions have similar range 
while at higher doping the charge correlation function 
has a much larger range compared to the spin correla- 
tion function. 

Two-dimensional case 

The phase diagram for the square lattice is shown in Fig. 
[5] The free energy for the square lattice as a function of 
the temperature is shown in Fig. [6] for U = 4, x = 0.1. 
The low temperature state is, as in the ID case, given 
by the phase x 7^ 0, A = (Red). There is a first order 
phase transition to the solution x / 0;A / (White) 
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U -2. U=4. U=5. 




Figure 5: Mean field phase diagram in the x — T plane as a function of U for a 16 x 16 square lattice. The colors represent the 
different types of solution minimizing the free energy: (1) x / 0, A = (Red); (2) x = 0, A / (Orange); (3) x / 0, A / 
(White); x> A = (Black). The calculated points are placed in the nodes of the mesh, the colors plotted between nodes are 
interpolated. 




Figure 6: Left Panel: Free energy for U = 4, x = 0.1 as a function of the temperature. The Black, White and Red dots follow 
the same color code as the phase diagrams. The Gray lines in background are the free energies of the rest of the converged 
solutions (most of which are unphysical, see main text). Central panel: Fermionic order parameters as a function of temperature 
for each phase separately, note that solutions that minimize the free energy are symmetric under rotations of the lattice. Right 
Panel: Minimum of the one particle excitation energies of the diagonalized bosonic and fermionic MF Hamiltonians. 



as the temperature increases. However contrary to the 
ID case there is only one second order phase transition 
for higher temperature where both hoppings and pairings 
vanish at the same time. Both phases (Red and White) 
present gapless bosonic (charge) excitations (see Fig. pi 
Right panel) but the minimum of the bosonic 1-particle 
excitations is located at k = {0, 0} and k = {tt, it} re- 
spectively. Furthermore the \ ^ 0, A = (Red) phase 
has no spin gap (see Fig. 6^Right panel, here we believe 
that the oscillations are due to finite size effects). Also 
note that contrarily to the one-dimensional case the fully 
incoherent phase has a fermionic (spin) gap smaller than 
the White phase. However, the White phase is the most 
stable one for dopings between 0.15 — 0.3 which indicates 
the presence of a spin gap (pseudogap) in this region. 

As a function of doping the free energy is minimal at 
half filling for the state x = 0;A/0 (Orange). This 
state is the zero doping limit of the state ^ / 0;A / 
(White) where the values of \ vanish. Moreover one 



observes a first order phase transition to the state x 7^ 
0, A = (Red) as the doping is increased. For U = 2 this 
arises away from half-filling resulting in two consecutive 
phase transitions with increasing values of x. the first 
from x = 0; A f (Orange) to X f 0; A ^ (White) is of 
second order, and the second from x ^ 0; A ^ (White) 
tox 7^ 0, A = (Red) is of first order. For U — 4, 5 only 
the first order transition is observed arising at x = 
between phases x = 0; A 7^ (Orange) and x 7^ 0, A = 
(Red). For still higher doping and sufficiently strong U 
there is another first order transition back to the state 
X^0;A^0 (White). 



C. Magnetic Properties 

By construction of the present MF decoupling there is 
no sub-lattice magnetization. Moreover since there is no 
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anti-ferromagnetic term in the Lagrangian the usual d- 
wave solutions of the t — J model are absent. Indeed, for 
two dimensions, the solutions with smallest free energy 
were found to be invariant under lattice rotation. How- 
ever, a collection of other solutions that break rotational 
symmetry of the lattice where found which are not the 
minima of the free energy, as shown in the Appendix. 
The Spin-Spin correlation (S. r .S r <) functions, at the MF 
level, are exponentially decaying with the distance and 
negative for r 7^ r' at finite T and for all non-trivial 
phases. This is compatible with the RVB picture where 
there is singlet formation for r and r' not nearest neigh- 
bors. 



D. Charge Properties 



Using that 



1 



Ck,, 



qSq+k.<T + <7,S q _ ff d q+k 



(20) 



we can write the Green's function as 

q 

{(4( T ) e l) \ S k+ q A' r ) s l+ q ,a) 

+ ( rf ^)<} (4-k,-Ar)s q -k,-,) Q 
+a (d q (r)e. q ) (4_ fc ,_»4- 9 , CT ) 

+(T ( e -q( T ) d q) ( S fc-g:<T( T ) S «-fe,-<r)o} 



(21) 



For the case with no pairing correlations between 
fcrmions and considering the bosons to be condensed 
leads to a Fermi liquid-like behavior for the electron 
Green's function where the pole like structure of the 
coherent part is determined by the spectrum of the 
fcrmions and the condensates give the Z renormalization 
factor. Since condensation was not allowed our phase 
diagram apparently misses two phases of Ref.[19j away 
from T = 0: the normal Fermi liquid phase and the su- 
perconducting (SC) phase. In one dimension the phase 
(x ^ 0, A = 0) has a charge (bosonic) gap that goes to 
zero with the temperature and a vanishing fermionic gap 
(see Fig. [3|. However, in two dimensions both gaps are 



zero in this phase (see Fig. [6]). Therefore at least at zero 
temperature there will be Bose condensation and assum- 
ing a small interlayer coupling we may have Bose conden- 
sation at finite temperature. Here we consider that the 
condensation temperature is smaller than the considered 
range of temperatures. Note that a conductor state is 
found for the x 7^ 0, A = phase as seen ahead in the 
spectral function where the spectral weight goes down to 
zero energy. For the other nontrivial phases a spin gap is 
present at low temperatures and one expects an insulator 
state. This is also corroborated by the spectral function 
analysis where no substantial spectral weight is seen near 
zero energy. 

No condensation implies that the "superfluid" density 
Pa(k) = V^ 1 \(a k dk} \ {a — d, e) is zero for all k and thus 
no SC correlations are present. Note that the symmetry 
of the SC gap is given by the SC correlation functions of 
the s particles. The superconductor correlation function, 
at the MF level, is given by 



,,2K 



K 



b k-KA b -k+K.-l J Q y 

x 44) 

\S-k+K,-l s k-K,l)o 



V 







(22) 



where 8k t o IS 27r periodic in k so the sum over K is reduced 
to two factors. However for the phases where A ^ the 
anomalous spin correlation functions (defined in Eq. 19 1 



are finite, which can be seen as a precursor for supercon- 
ductivity. 



IV. SPECTRAL FUNCTION 

An efficient way to obtain information about the ex- 
citation spectrum of a strongly correlated system is 
through the spectral function. This is defined as the 
imaginary part of the Green's function and is directly 
measurable through photoemission experiments. For the 
theory considered here the electron spectral function at 
the MF level can be written as 



A a {u>,k) 



^ E [ S ( UJ + W e, 9 + U)F,-q-k) [1 + nb,e{l) ~ n f(~ ( l ~ \ U F,q+kVB ,-q + U B ,-qVF,q+k\ 

q 

+6 (oj - oj d . q + 0J F , q -k) [n/(q - k) + n b4 (q)] \u B . q u F . k _ q + v F ^ q v B ^ q \ 2 
+5 (oj - oj d . q - 0J F ,k- q ) [1 + n b4 {q) - n f (k - q)] \u Fk - q v B ^ q - u B ^ q v F ^ k _ q \ 2 
+6 (u + u! etq - u F>q+k ) [n f (q + k) + rib,e(<l)] \uB,- g u F>g+ k ~ v FA+k v B .- q \ 2 



(23) 
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Figure 7: Spectral function and single particle energies for the Hubbatd chain: (a) Phase \ 7^ 0; A = computed for 
U = 4; x = 0.3; T = 0.06. (b) Phase x / 0; A / computed for U = 4; x = 0.5; T = 0.06. (c) Phase x = 0; A / computed 
for (7 = 4; x = 0; T = 0.06. The U — x — T parameters were chosen such that each solution is a minimum of the free energy. 
Colored lines are plotted along regions of considerable spectral weight and are obtained as follows: (Dark Blue, dotted line) 
varying the momentum of the fermions while keeping the momentum of the bosons to its minimal energy; (Green, dashed line) 
fixing the fermion energy to its minimal value and changing the momentum of the bosons; (Red, solid line) requiring that both 
excitations have equal velocities. 
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and has information about the excitation spectra and 
their spectral weights. 

Experiments for one-dimensional insulators and con- 
ductors have shown the fractionalization of the electronic 
degrees of freedom inside the strongly correlated system. 
For instance experimental results for the one-dimensional 
conductor TTF-TCNQ have been interpreted using the 
Bethe ansatz solution showing clearly the traits of the 
spinon and holon branches associated with the spin and 
charge degrees of freedom (see Fig. 9 in [S] and Figs. 1 
in (9j HI] )• In the context of the high-temperature su- 
perconductors similar results have been obtained for the 
spectral function in |31l I3"2] . These results have also been 
interpreted in terms of some fractionalization of the de- 
grees of freedom in a way similar to the one-dimensional 
case. 

Fig. [7] shows the spectral function for the one- 
dimensional case computed at the MF level for the three 
non-trivial MF solutions. The regions of larger spectral 
weight correspond to three different processes. Two of 
those processes can be interpreted as the fractionaliza- 
tion of the electron since the fermions or the bosons (Blue 
or Green respectively) are fixed to the minimum value of 
their energy and the other species is allowed to move 
along its energy band. The high spectral weight of such 
regions is also observed in the spectral function computed 
based on the exact solution of the Hubbard model in ID 
(REF). These are the so-called spectral lines "c" and "s" 
in refs. [3 [HI \H\- The third process can be interpreted 
as propagation of an "electron-like" degree of freedom 
(Red) since it corresponds to a boson and a fermion with 
the same velocity propagating together. These processes 
typically define the boundary lines of the region with a 
(nearly) non- vanishing spectral weight. 

Consider first the top panel of Fig. [7] and the negative 
energy region corresponding to photoemission. In the re- 
gions of higher spectral weight at low energies we see that 
the lowest energy branch (s branch) is obtained fixing the 
bosons at their lowest energy and changing the fermionic 
(spin) particles along their bands. It is therefore a spin 
branch. Below this line at higher energy (recall that by 
definition u> < for photoemission and therefore higher 
energies are more negative) there is a line that is obtained 
fixing the fermion at the Fermi surface and changing the 
boson energy along its band. It is therefore a charge 
(holon) band (c branch) . These two lines are clearly sep- 
arated as in the exact description from the Bethe ansatz 
and the experimental results [TT]. They merge at the 
Fermi level as in the exact solution. Note that in the re- 
gion of high spectral weight there is a contribution from 
a line obtained taking the velocities of the fermions and 
bosons as equal. This means an "electronic-like" excita- 
tion. For larger momenta spin and a charge branch also 
emerge as in the exact solution. We also show the in- 
verse photoemission spectra (positive energies) including 
the lowest and the upper Hubbard bands. In this phase 
there is a finite spectral weight at the Fermi energy which 
implies a conducting phase. This is consistent with the 



gapless fermionic band and the nearly gapless bosonic 
lowest band. 

In the middle panel we consider the phase where both 
the hoppings and gap functions are finite (White phase). 
In this case the spins have a gap and the bosons are nearly 
gapless. The spectral function has now a pseudogap at 
the Fermi level since there is a small spectral weight. The 
finite energy structure is however quite similar to the fully 
conducting phase. Note that another contribution to this 
region is obtained fixing the fermion energy to its minimal 
(finite) value and changing the momentum of the bosons 
(nearly gapless). Note again that in the region of high 
spectral weight there is again a contribution from a line 
obtained taking the velocities of the fermions and bosons 
as equal. 

Finally, in the lower panel we consider the Orange so- 
lution where the hopping parameters vanish and the gap 
functions are finite. In this case both the fermions and 
the bosons are gapped and the spectral function has a 
large gap. Also, the results are presented at half-filling 
and therefore the system is an insulator, as expected from 
the exact solution. 

The spectral function and density of states 

A(w) = ^~£4,(w,fc) (24) 

k,(T 

for the square lattice are shown in Figs(8]and[9] For the 
X 7^ 0, A = (Red) phase it presents the same quali- 
tative features as the one dimensional case. Namely a 
strong spectral weight is observed due to excitations cor- 
responding to an empty site at k = and a spinon that 
carries the momentum of the physical electron (s branch) . 
Excitations for which the spinon is taken at the Fermi 
surface and the holon carries the difference of momenta 
presents also some spectral weight (Green lines of Fig|9j . 
Figj9] shows Green lines (c branches) corresponding to 
the spinon momentum at the Fermi surface in the direc- 
tion (0, 7r), (tt, tt), (tt, 0) directions. The line reaching the 
Fermi energy corresponds to a fermion with momentum 
in the (0, 0) — (tt, tt) segment. This scenario of two spec- 
tral lines leaving the Fermi-surface of the spinons is also 
obtained by other recent approaches [22]. The results for 
the White phase show a pseudogap structure as shown 
in the density of states A(u>). 



V. DISCUSSION 

We have extensively explored translationally invariant 
MF solutions of the Hubbard model for a one dimensional 
chain and for a square lattice using an electron represen- 
tation introduced in |13| and a MF decoupling in terms 
of link variables. In two dimensions this includes non- 
trivial symmetries of the MF solutions as well as nematic 
(translationally invariant but not rotationally invariant) 
and quasi-lD phases (but no flux phases). Despite all the 
freedom in the choice of solutions, the ones that minimize 
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Figure 8: Spectral function and density of states computed at the MF level for a square lattice 16 x 16. Upper panel: phase 
X 0, A = (U = 4,T = 0.01, x = 0.1). Lower panel: phase x + 0, A / (U = 4, T = 0.2, x = 0.1) 
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Figure 9: Spectral function for the square lattice for the phase 
X 7^ 0, A = along the nodal line. The lines green are ob- 
tained fixing the fermion momentum at the Fermi level along 
the directions (1) = (tt,0); (2) = (7r/2,7r/2); (3) = (0,tt) and 
varying the momentum of the bosons. The blue lines are ob- 
tained fixing the bosons energy to be minimal (corresponding 
to zero momentum) and varying the fermionic momentum. 



the free energy were found to be invariant under lattice 
rotations. However, in some phase space regions, the free 
energy difference between these less conventional phases 



(non rotationally invariant) and the symmetric ones (ro- 
tationally invariant) where found to be quite small sig- 
naling a possible stabilization of such phases by some ex- 
tra coupling in the Hamiltonian. The symmetric phases 
where classified according to the MF order parameters 
and their physical properties where obtained. 

Generically we found a gapped phase for both 
fermionic and bosonic (spinons and holons) degrees of 
freedom for zero doping. Contrarily to the 2D case, where 
this phase was only found for half filling, in the ID case 
this phase extends to finite doping at finite temperatures 
and its size in the T — x phase diagram decreases with 
U . This is in no contradiction with the exact solution for 
the ground state since for zero temperature only at half 
filling have we found this phase; the system is a conduc- 
tor away from half filling. Although no magnetic order is 
obtained from the present MF ansatz this state is clearly 
a Mott insulator as one expects at half filling. 

A conducting phase was found to be dominant for small 
and moderate doping from zero to quite high tempera- 
tures presenting the qualitative features of a RVB state. 
For the one dimensional conducting phase one can iden- 
tify some of the features of the spectral function of Fig. [7] 
with the ones obtained from the exact Bethe Ansatz so- 
lution. In particular the lines carrying the most part of 
the spectral weight can be identified with three kinds of 
excitations. Two of them corresponding to a creation or 
annihilation of slave particles with different velocities. It 
is tempting to interpret such lines with fractionalization 
of the initial degrees of freedom. In the two dimensional 
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Figure 10: This figure shows how the results displayed in Fig. [TT] should be interpreted. 



case these lines are also obtained although they are not 
so clearly defined, in the sense that the spectral weight is 
not very pronounced for the c branches. Near the Fermi- 
energy such two lines are compatible with the description 
in [22 . Lines where both slave particles have equal ve- 
locities are possible to draw only in the one dimensional 
case since they correspond to a surface in 2D. In this 
case they typically represent boundaries for the spectral 
weight. As in the former case one can try to interpret 
these lines as "electron- like" particles since they repre- 
sent states where a slave boson and slave fermion travel 
together with the same velocity. In the region where the 
X / 0,A / phase appears the temperature is larger 
than the spin gap, the charge gap being zero. We note 
that in the one-dimensional case in the conducting phases 
the fermions are gapless and the bosons are gapped while 
in the square lattice the fermions are gapped and the 
bosons are gapless. 

The results presented here describe the finite-energy 
finite-temperature phase diagram. At very low tempera- 
ture a low energy theory where the lowest energy bosonic 
branch is condensed and the higher energy bosonic 
branch is frozen, will be presented elsewhere. It is par- 
ticularly interesting in this low energy regime to consider 
frustrated lattices where it is expected that fractional- 
ization will appear for either a conducting system, or for 
insulating systems where the possibility of spin liquids 
has been proposed. 
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Appendix A: Mean Field Solutions 

In this section we present a table with some of the 
solutions found for the square lattice, solving the MF 
equations. These solutions have low values of the free 
energy even though not absolute minima. 



The results presented in Fig. 11 are obtained for U = 
4, x = 0.2,T=0.1 for the two dimensional square lattice 
and are representative of the types of solutions found 
in other points of phase space. Solutions (1),(3) and 
(4) are the ones considered in the main text since they 
alternatively minimize the free energy depending on the 
considered phase space region: (1) corresponds to the 
X 7^ 0, A — phase, (3) corresponds to the x^0,A^O 
phase and (4) is the incoherent phase where \ — 0, A = 
0. Note that the solution (2) referred in the main text 
corresponds to the zero doping limit of solution (3), the 
one shown in Fig. [TTJ is also a solution of the type \ — 
0,A ^ but it never presents the lowest free energy. 
The table entries are complex plane plots representing 
the MF parameters: columns 2 to 5 represent the values 
of Xf/b.8i Af/b^ for 8 — e x (orange arrow) and 5 = e y 
(black arrow); column 6 displays the values of (i (black 
arrow) and iX (orange arrow). The scales of the axes are 
shown in Fig. [10] We note some remarkable solutions 
obtained here: solution (4) is an anisotropic \ 7^ 0, A = 
solution; solutions (9) and (10) are related by a lattice 
rotation and represent one dimensional- like correlations 
where in one of the directions the MF parameters are 
zero; solution (16) has x 7^ 0, A = in the x direction 
and x — 0, A ^ in the y direction. 
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